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Abstract We give an expression for the Lindblad torque acting on a low-mass planet 
embedded in a protoplanetary disk that is valid even at locations where the surface 
density or temperature profile cannot be approximated by a power law, such as an 
opacity transition. At such locations, the Lindblad torque is known to suffer strong 
deviation from its standard value, with potentially important implications for type I 
migration, but the full treatment of the tidal interaction is cumbersome and not well 
suited to models of planetary population synthesis. The expression that we propose 
retains the simplicity of the standard Lindblad torque formula and gives results that 
accurately reproduce those of numerical simulations, even at locations where the disk 
Qh. temperature undergoes abrupt changes. Our study is conducted by means of cus- 

tomized numerical simulations in the low-mass regime, in locally isothermal disks, 
and compared to linear torque estimates obtained by summing fully analytic torque 
estimates at each Lindblad resonance. The functional dependence of our modified 
Lindblad torque expression is suggested by an estimate of the shift of the Lindblad 
resonances that mostly contribute to the torque, in a disk with sharp gradients of 
temperature or surface density, while the numerical coefficients of the new terms are 
adjusted to seek agreement with numerics. As side results, we find that the vortensity 
related corotation torque undergoes a boost at an opacity transition that can counteract 
migration, and we find evidence from numerical simulations that the linear corotation 
fNj , torque has a non-negligible dependency upon the temperature gradient, in a locally 

' isothermal disk. 
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1 Introduction 

As the number of extrasolar planets discovered increases continuously, there is a 
growing volume of models of planet population synthesis, in which one evolves a col- 
lection of planets in a model of protoplanetary disk, s o as to compare the statistics of 
real pla n etary systems to thos e of the model outcomes dLvra et allbOlOtlMordasini et all 



2009allbl Ida and Lin, 2008k a). One key feature of models of planetary population 



synthesis is the tidal torque that exists between a planet and its parent protoplanetary 
disk, and which drives a radial migration of the planet. Since this torque is large, plan- 
ets migrate a substantial fraction of their initial semi major axis over the disk lifetime, 
and a model's outcome strongly depends on the underlying torque expression. The 
tidal torque has two components: 

- The differential Lindblad torque, exerted by a wake excited by the planet in the 
neighboring flow, where the latter becomes supersonic with respect to the planet, 
owing to the strong Keplerian shear. 

- The corotation torque, which corresponds to an exchange of angular momentum 
between the planet and the coorbital material. 

Early efforts have intensively focused on the differential Lindblad torque, which is 
no w known in great detail in power law disks. The well known standard formula 
of Tanaka etal (2002) has proved to account for the torque measured in numeri- 



cal simulations with remarkable accuracy. More recently, the corotation torque has 
received a detailed attention, and was shown to have two components: one related 
to the vortensitjQ gradie nt across the orbit, and ano t her one related to the entropy 



i 1 1 1 1 

gradient across the orbit (Baruteau and Masset, 2008; Paardekoorjer and Parjaloizou, 



2008 



Casoli and Masset, 2009; Masset and Casoli, 2009, 2010; Paardekooper et al, 



2010). The corotation torque, moreove r, depends on the dissipative processes at work 



in the disk. Masset and Casoli (2010) have provided general corotation torque for- 



mulae that may be applied to planets embedded in disks with arbitrary vortensity and 
entropy profiles, and with arbitrary values of the viscosity and thermal diffusivity. 

The domain of validity of the formulae for the corotation torque now supersedes 
that of the Lindblad torque: the corotation torque formulae apply to arbitrary profiles, 
and only requires the evaluation of radial derivatives at the orbit, whereas the Lind- 
blad torque standard formula exclusively applies to disks where the surface density 
and temperature profiles are power laws of the radius. When this is not the case, 
one has to resort to basics and evaluate the tor que individually at each Lin dblad 
resonance. This was the approach undertaken by Menou and Goodman (2004), who 
stressed the high sensitivity of the torque to the resonance location. Their key result 
was that the Lindblad torque can be substantially reduced at locations, in the disk, 
where the local profiles of surface density and/or temperature cannot be approxi- 
mated by power laws, namely at opacity transitions. 

This approach, nonetheless, would be very costly in planetary population synthe- 
sis models, as it would require the evaluation of a series at each time step, for each 
planet. To the best of our knowledge, all the planetary population synthesis models 



The vortensity is here defined as the ratio of the vertical component of the flow vorticity over the 
disk's surface density. 
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performed thus far have made use, in different flavors, of standard torque formulae for 
type I migrating objects, in some case including a uniform, ad hoc reduction factor 
so as to moderate the drastic effects of standard type I migration. 

The present work is an attempt to provide a simple generalization of the standard 
Lindblad torque formula, that can be applied in arbitrary disk models (/. e. where 
the profiles of surface density and temperature are not power laws of radius). The 
reader may glimpse at Eq. ( 179} , which constitutes our main result, and compare it 
with Eq. (TTOl l. which is the standard expression in power law disks. As can be seen, 
there is very little additional complexity in evaluating Eq. d79l l, which can be readily 
implemented in models of planetary population synthesis. 

We firstly define the notation in section|2] We then discuss in some detail the rele- 
vance of numerical simulations in the low-mass regime in assessing the linear torques 
in section |3] and we present standard torque expressions in section [4] In particular, 
we need to calibrate the corotation torque in order to properly subtract it from our 
runs, so that we get exclusively the Lindblad torque. We present our method in sec- 
tion |5] We then derive the effects of higher order derivatives of the surface density 
or temperature on the Lindblad torque in section |6] As we shall see, the end result 
is particularly simple, as it turns out that only the third derivative of the temperature 
has a significant impact on the Lindblad torque. We then apply this formula to a disk 
model aimed at representing an opacity transition in section [7] where we present the 
generalized torque expression in definitive form. In that section we also estimate the 
corotation torque. We find that at the transition the corotation torque also features a 
term proportional to the third derivative of the temperature, with a sign opposite to 
that of the Lindblad torque. We discuss our results in section|8]and draw our conclu- 
sions in section|9] 



2 Notation 

We consider a planet of mass M p orbiting a star of mass M» on a circular orbit of 
radius a, with orbital frequency Q p . The planet is embedded in a gaseous disk with 
surface density S{r) and temperature T(r). The vertical scale height H(r) of the disk 
is related to the gas sound speed c s (r) and to the local Keplerian frequency Qk(/) by: 

c s (r) 

We oftentimes use the disk aspect ratio or relative thickness defined as h(r) - H(r)/r, 
and the planet to star mass ratio q = M p /M*. We denote to the vertical component of 
the flow's vorticity: 

CO = -d r (r 2 Q). (2) 
r 

We define hereafter notation specific to non power law disks. We denote with 
a quantity related to the i th derivative of the surface density: 

a| = ( _ ir NH)AW^gft (3) 
(dlogr)' 
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We thus have: 



a] 



<*2 



d\ogZ 
d log r ' 
d 2 logr 
dlogr 2 ' 
</ 3 log 27 



(4) 



(5) 



(6) 



d log r 5 

In a power law disk, all with i > 1 cancel out and we have Eq cc r _ct| . Our 
value of a\ coincides with the standard notation a of iTanaka et all d2002l) . The fac- 
tor (-l) mm ('' 2 ) in Eq. (0 is exclusively meant to conserve the minus sign for the first 
derivative, so as to save the standard practice. In contrast, higher order derivatives do 
not feature this minus sign. In a similar fashion we denote with /?, a quantity related 
to the i derivative of the temperature: 



= (-l)™ n C'. 2 V' 



d' log T 



(7) 



(dlogr)' 

Our B\ therefore coinci des with the value of C o f iD'Angelo and Lubow (1201 Oh . or 
with the value of B of iPaardekooper et al d2010l) . but it is twice the value of f5 of 
Tanaka et all d2002l) . where /3 is the slope of the sound speed rather than the temper- 
ature. We work out how the differential Lindblad torque behaves in a disk in which 
some of the a, and/or (z > 1) do not cancel out. We limit ourselves to i < 3. As 
we shall see in section [TJ there is no need to push the expansion further. We entertain 
separately the role of the a, and (z > 1), assuming that the resonance shift is suffi- 
ciently small that the global effect is obtained by adding up individual contributions 
of higher derivatives of surface density and temperature. We shall often refer to the 
residue of the surface density profile with respect to the power law surface density 
with same a\ as a perturbation (the same holds for perturbations of temperature). 
It should be noted however that this perturbation is not related to the planet. It is 
axisymmetric, and the "perturbed" disk is in rotational equilibrium. 



3 Relevance of numerical simulations in assessing the torques acting on a 
low-mass planet 

Since full hydrodynamical simulations are non-linear in nature, it is questionable 
whether the torques obtained in simulations in the low-mass regime can be regarded 
as estimates of the linear torque value. As a prerequisite, it is therefore necessary 
to confront the torque thus obtain ed to the torque obtained independ ently by a lin- 
ear analysis. This is the work of IPaardekooper a nd Papaloizou (2009). Part of their 



findings (the ones that we shall need in this work) are as follows: 

1 . For low-mass planets, non-linear effects are found at all masses, and identified 
with the onset of the horseshoe drag regime, provided the disk's viscosity is suffi- 
ciently low. In that case, the linear corotation torque exists only temporarily upon 
insertion of the planet. The timescale for the onset of non-linear effects scales 
with q~ l/2 . 
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2. In the inviscid case, the (non-linear) corotation torque saturates at larger time, 
leaving only the linear Lindblad torque. 

3. At early time or at any time for large viscosity, the total torque acting on the planet 
is found to be in excellent agreement with the total linear torque. 



Paardekooper and Papaloizou (2009) established these results in disks with a locally 
isothermal equation of state, that is disks in which the sound speed or temperature is 
an arbitrary function of r, constant in time. Admittedly, they obtained these results by 
considering various temperature and surface density profiles which are exclusively 
power laws of the radius. 

We shall use the result 2 above in the very same framework (see appendix lAl). 
i. e. for power law disks, in order to establish an expression for the Lindblad torque 
in these disks, that is valid for our particular setup (two dimensional disks and a 
planetary potential softened over a given length.) In regard to the above we shall 
consider that expression as a good approximation of the linear Lindblad torque. 

We shall use the result 3 also in power law disks (see appendix|A]i, so as to get an 
expression for the total torque, which, again, we regard as a good approximation of 
the total linear torque for the reasons listed above. This enables us to get an expres- 
sion of the corotation torque at early times (by subtracting the Lindblad torque), as a 
function of the gradients of vortensity and temperature. This step is needed to enable 
us to subtract the corotation torque from the runs in which the disk profiles are not 
power laws of radius, so as to get the Lindblad torque. 

Finally, we shall measure the total torque at early time in calculations where the 
surface density or temperature profiles are not power laws of the radius. Although, 



strictly speaking, this is not contemplated in the framework of Paardekooper and Papaloizou 
d2009h . it is reasonable to consider that the torque acting on the planet has a linear 
value at this stage. Regardless of the profiles, the response of a disk to a planet of 
low mass {q/h 3 <K 1) should be linear over a time interval that can be made arbitrar- 
ily large by decreasing the planetary mass. We further comment that in our case, we 
have q/h 3 = 1 .6 ■ 10~ 2 , almost an order of magnitud e smaller than that of the fiducial 



calculation of Paardekooper and Papaloizou (2009), which is q/h 3 =0.1 



We conclude that numerical simulations, albeit non-linear in essence, can be used 
to assess with a good accuracy the linear torques acting on a low-mass planet. The 
de gree of accuracy can be e stimated either by comparing our data to the expression 
of Paardek ooper et all (120101) for power law profiles, or by evaluating the fitting error 
of our data (see appendix|A]l. From both methods one can infer a degree of accuracy 
of a few percents, which is far enough for our purpose. 



4 Standard torque expressions 

Standard Lindblad torque formulae, be they two or three dimensional, all share the 
same fornj|: 

r L = -(k +k 1 a 1 +kiJ3 1 )r Te{ , (8) 



- As a convention throughout this work, all torques refer to the torque exerted by the disk onto the 
planet. 
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where kg, k\ and k^ are numerical constants of order unity, and where r re f is the torque 
normalization factor given by: 

T ref = Za 4 Q 2 p q 2 h- 2 . (9) 

We refer to Lindblad torque formulae having the form of Eq. (|8]) as standard to ex- 
press that they exclusively involve a linear combination of the surface density and 
temperature first derivatives at the position of the orbit. Strictly speaking, standard 
Lindblad torque formulae are valid only in disks in which both surface density and 
temperature are power laws of the radius. 



The expression of the Lindblad torques alone are not explicitly given by Tanakaeta] 



(2002), but they can be inferred from their Tab. 1 and 2 (see iMassetL 20081) . This 
yields {ko,k\) = (3.2,1.468) in the two dimensional case without softening, and 
(ko,ki) = (2.34, -0.099) in the three dimensional case. Note that ki was left unex- 
plored in this original work, as the authors restricted themselves to globally isother- 
mal disks (i.e. without a temperature gradient). 

In the present work we shall deal with two dimensional disks and a softened 
planetary potential (we chose a softening length of the potential e = 0.6H). For our 
purpose, we therefore have to firstly establish the standard Lindblad torque formulae 
in those disks. This is the purpose of appendix [A] in which we evaluate simultane- 
ously the standard torque expressions for the Lindblad torque and for the total torque. 
These expressions are respectively 



and 



r num,l = _ (2 QQ _ Q l6ai + inp^f^ (10 ) 



r t n u t mJ = -(1.09 + 0.45ai + 0.5300^. (11) 



where the superscript (which stands for numerical, linear) conveys the method by 
which the expressions have been obtained, that is by means of numerical calcula- 
tions, as an approximation to the linear value, along the lines of what has been ex- 
posed in section [3] We n ote that Eq. ( TTOb is in good agreement with Eq. (14) of 



Paardekooper et all (12.010b . It is of interest to confront Eq. (TTOb to a fully analytic 



torque estimate (see appendix|B]for details), which yields: 

r Unear _ _ ( g _ q + 2.09/?,)/^. (12) 

This estimate differs significantly from Eq. ([Tot , which is also very close to the esti- 



mate obtained by solving numerically the linear equations (IPaardekooper and Papaloizou 



2009). Similar discrepancy between fully an alytic results and resu l ts obt ained by 



other methods have already been noticed by iD'Angelo and Lubow d2010h . Still, it 



broadly captures some of the salient features of the torque: it is negative for reason- 
able values of the surface density or temperature gradient, a faster decay outwards 
of the temperature (larger f}\) shifts the resonances inwards and renders the torque 
more negative (i.e., the coefficient of fi\ is negative), and for similar surface density 
gradients and temperature gradients (which shift the resonances in almost exactly 
the same amount), there is, in the case of a surface density gradient, an additional 
weighting by the local surface density which tends to counteract the resonance shift 
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(i.e., the coefficient of a\ is larger than that of f3\). However, whereas Eq. (fTUt shows 
that both effects nearly cancel each other, as th e coefficient of a\ almost vanishes [an 
effect known as the pressure buffer, shown bv iWardl (119971) 1. it is not quite the case 
of Eq. ( fl2l i. A near cancellation is found, in fully analytic torque estimates, only for 
very low values of the smoothing length, and the coe fficients of a\ and P\ are foun d 



to depend heavily on the smoothing length (see also Me nou and Goodm an, 2004). 



The discrepancy between a fully analytic estimate and the true linear torque value 
(given either by the numerical solution of the linear equations or estimated by means 
of adequate numerical simulations) likely arises from the approximation of the torque 
cut-off in fully analytical estimates. This cut-off is given by the Eq. d9lT ), which, to- 
gether with the expression for r m , the location of the resonance, determines how the 
Lindblad torque behaves for highe r order resonances. The denominator of Eq. ( l9TT l, 
in particular, is only approximate ( Artvmowiczl 1993 ). 



The expression of the corotation torque, given by: 

j^num.l wium.l wium,l /rv m /: 1 , n c o a \ / 1 o \ 

r c = r tot r L = (0.91 - .61ai + 0.58j8i)r ref (13) 

must be recast for our needs to: 

r™ mJ = (0.61<V + 0.58/3, )r ref (14) 

where *V = dlog(X/oj)/dlogr is referred to for brevity as the vortensity gradient. 
Eq. ( fT3l ) and (fl"4b are equivalent in power law disks. However, it is the form of Eq. (fT4l 
that must be used for arbitrary profiles, as the vortensity gradient *V can differ sig- 
nificantly from 3/2 - a\, as we shall see in section 1741 We shall come back to this 
substitution in section Rx2l where we compare the effect of higher order derivatives on 
the torque obtained either by a fully analytic approach or by numerical simulations. 
We shall see that these effects, which affect preferentially lower order resonances, are 
not plagued by the inaccuracies due to the torque cut-off. This allows a direct com- 
parison of analytics and simulations, which give essentially the same results, thereby 
validating Eq. d!4l . 

It is instructive to note that the coefficie nts of a\ and J3\ in E q. ([T4l are approxi- 
mately equal. One can notice that Eq. (11) of lTanaka et all d2002h is suggestive of this 
result. Getting rid of all terms of this equation that involve derivatives in z, which do 
not exist in a two-dimensional disk, one sees that the pole at corotation for rj' m consists 
of two contributions: one in the vortensity gradient, and another one that scales with 
2/3 (remember that Tanaka et al's ft is half our/3i). Eq. (fT4l will be used from now 
on to subtract the corotation torque (again considered here as a good approximation 
to the linear corotation torque), from the total torque at early times in calculations 
where the disk profiles are not power laws, so as to get the Lindblad torque. 



5 Assumptions and set up 

5.1 Wake shift 



In disks in which the surface density profile or the temperature profile is not a power 
law of the radius, the locations of Lindblad resonances differ from their locations in 
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a power law disk with same value of a\ and The Lindblad torques, which sen- 
sitively depend on the resonance location, can therefore undergo large changes, not 
captured by a standard torque expression. Our approach is threefold: we give an an- 
alytical order of magnitude of the effect, we check our predictions by means of ded- 
icated numerical simulations, and by means of a completely analytic estimate given 
by summing the contribution of all Lindblad resonances (see appendix[B]) Although 
this last method yields results more accurate than our order of magnitude estimate, 
it does not provide much insight into the relationship between the local properties of 
the density and temperature profiles on one hand, and the Lindblad torque properties 
on the other hand. For this reason we find it relevant to get an estimate of magnitude 
of the effects of higher order terms as follows: we focus on the resonances that most 
contribute to the torque, which correspond to azimuthal wave numbers m = OilT 1 ). 
For such values of m, we evaluate the resonance shift, which we improperly call the 
wake shift. In order to translate this wake shift into a torque variation, we firstly cal- 
ibrate the one-sided torque dependence on the wake shift by a series of dedicated 
calculations, in which we evaluate how the (one-sided) Lindblad torque changes as a 
resonance is shifted, by tuning the orbital frequency of the planet that remains on a 
circular orbit of fixed radius. This is exposed in appendixO From the standard torque 
formula at a given resonance, given at Eq. (I90K we see that under a perturbation of the 
disk profiles the torque can change either as a result of the forcing potential (hence a 
shift, in our case) or as a result of a change of the "width" of the resonance, given by 
the denominator in rdD/dr. The latter effect is negligible in our case, and most of the 
torque change is accounted for by the resonance shift. 

We note that our simplified method of evaluating the "wake shift" presents an 
obvious flaw: the major contribution to the Lindblad torque change may not arise 
from the resonances that most contribute to the Lindblad torque. This method gives 
nevertheless a reasonable order of magnitude of the torque variation, and, more im- 
portantly, it suggests a functional form for the torque variation due to higher order 
derivatives. Once this form has been obtained, one can resort to numerical simula- 
tions to adjust numerical coefficients. 

We comment that even though we can obtain the linear Lindblad torque by a lin- 
ear analysis, which is self-contained and does not feature free parameters, the result is 
not strictly amenable to a simple form such as the one we are aiming to. An obvious 
reason for that is that the higher order derivatives of the di sk profiles at the location of 
each Lindblad resonance feature in the result. The reduction of the Lindblad torque 
to a simple approximate expression still involves the identification of a physically 
grounded, simple functional form, and the subsequent adjustment of free parame- 
ters. Our final expression should therefore be regarded as an empirical estimate, even 
though its functional dependence is motivated by considerations on the wake shift. 



5.2 Set up and method for numerical checks 

We undertake in the following two kinds of numerical experiments: 
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1 . in sections 16.21 and 16.31 we check the dependence of the Lindblad torque on a 
given higher order derivative of the temperature (section 16.21 ) or surface density 
(section l6.31 l. 

2. In section|7]we examine the Lindblad torque behavior at an opacity transition. 
For the first series of numerical experiments, we adopt disk profiles that are power 
laws of radius, except on the annulus a — E<r<a + E. Over this annulus, the 
value of f3i or a,- of interest is constant. The disk profiles outside of the annulus are 
extrapolated so that a\ and /3\ are continuous over the whole computational domain. 
We chose E = 2H or E = 3H, in order to cover the domain over which the exchange 
of angular momentum between the disk and the planet takes place. 

All numerical simulations were performed with the FARGO code0 ( M asset , 2000alfb" ), 
in its locally isothermal version. In all cases we evaluate the linear total torque by a 
time average from 2 to 4 orbital periods, then we subtract the linear corotation torque 
expression given by Eq. (fl4l . An implicit assumption underlying this procedure is 
that the linear corotation torque does not depend upon higher order derivatives of the 
temperature or vortensity. This is a reasonable assumption, as this torque is associ- 
ated to a singularity at corotation, whereas the Lindblad torque owes its sensitivity to 
higher order derivatives to its non-local character. Note however that the corotation 
torque must be cast in terms of the gradient of vortensity and gradient of temperature. 
The former, in particular, can bear a dependency on higher order derivatives of the 
surface density or temperature. 

The calculations of sections 16.21 and 16.31 were performed with a planet of mass 
M p = lO^M,, embedded in a disk with aspect ratio h = 0.05 at the planetary orbit. 
The mesh extends from/? m i n = 0.6 to 7? max = 1.4, and has resolution Af ra d = 1200 and 
Ng - 1500. Non-reflecting boundary conditions were used, but those should have no 
impact on the results, since during the early stage in which we evaluate the torque 
the wake has hardly reached the mesh boundaries. More details about the setup of 
section Q will be given therein. 



6 Individual effects of higher order derivatives 

6.1 Expression of the wake shift 

The "effective" location r m (e.g. lArtvmowicz , Il993h of a Lindblad resonance asso- 
ciated to the component of the planetary potential with m-fold symmetry is given 
by: 

D(r m ) = 0, (15) 

where 

2 2 
III ( 

D(r) = k 2 - m 2 (Q - Q„f + — (16) 

r l 

The effect of a perturbation of the sound speed or surface density profile is to induce 
a variation 6D of D, which shifts the resonance location by the amount: 

5D 

5rm = 'T^TTT- (17) 
(dD/dr),. m 
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The rotational equilibrium of the disk reads 



I 

while the expression of the epicyclic frequency is 



rti 2 = rQ\ + — , (18) 



K 2 =AQ\Q+-r^- . (19) 



\ 2 dr 



Using Eqs. ( fT~8b and ( fT9l we infer: 



2 ^ A8r P ^ d i d >P\ ,->n, 



r i7 dr\rS 

which can be recast, in a power law disk, as: 

^ 2 = ^[l+(a l + J 6 l )( y 6,-2)/ 1 2 ]. (21) 

Deriving D with respect to r, we get: 

<9D 3£ 2 3Qlm . „Ql 

^- = + — 5- - (2 +J 6,)m 2 / ; 2 ^, (22) 

or r r r 

where henceforth the upper (lower) sign is for the outer (inner) Lindblad resonance, 
and where we have retained terms in m 2 h 2 but not those in h 2 . Since the components 
that most contribute to the Lindblad torque have m ~ hr l » 1, we can simplify, for 
our purpose, Eq. (l22l as follows: 

3D 3Qlm 

-a- = (23) 
or r 

from which we infer 

rSD 

6r m = ±—3— • (24) 

We contemplate below two distinct cases: (i) the case for which the rotation profile is 
altered by a perturbation of temperature, the surface density remaining constant, and 
(ii) the opposite case, in which the profile of temperature remains unchanged, but a 
perturbation of surface density is introduced. 



6.2 Case of higher order temperature derivatives 

The variation of D under an axisymmetric perturbation of sound speed (or tempera- 
ture) is: 

6D = jd r 6c 2 + rd r (d r Sc 2 ) - 2m 2 (Q - Q p )6Q + '-^dc 2 = 0, (25) 

where 5Q is the perturbation of angular velocity which can be related to the pressure 
perturbation by means of Eq. (fT8l : 

d r 5c 2 

dQ =iaF- (26) 
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Using Eq. d26b . Eq. (f25b can be recast as: 



„, , md r 6c 2 m 2 
6D = d\Sc 2 ± — — + -r6c 2 , 
r r l 



(27) 



where we retain only the highest order terms in m. The radial shift of the resonance 
is thus: 



Sr„. 



2>Ql \ m ' 



6c, + d r 6c 2 



m 9 \ 



(28) 



We now specialize to the case of a non-vanishing fa. Exploiting the fact that fa 
is non-negligible whenever fa varies over a length scale much shorter than r, we can 
transform Eq. (0 into: 

hr 2 d 2 6c 2 



fa^ 



dr 2 



Denoting with x„ - r m - r c , where r c is the corotation radius, we can write: 



Sr m = 



c 2 fa 



3Q\hr 2 \ m 



r mx m 
+ — +x m + — — 
r 2 



(29) 



(30) 



We denote m max the value of the azimuthal wavenumber for which the torque is max- 
imal. Since m max is of order ~ 2h~ l /3, we substitute the location of the resonance x m 
with +AH, where we introduce A a 1 as a free parameter. This yields: 



Sx H 



fah 
3 



(uA) 1 + 1 + -uA 
2 



(31) 



where u = m m . dx h, and where x w = r w - a is the distance of the wake to the orbit. 

We can fur ther transform this expression using the value of u inferred for instance 
from Fig. 5 of Masset (2008), which gives u a 0.63. Similarly we take A = (2/3)(l + 
u- 2 ) l/2 « 1.25 dArtvmowiczlll993l) . This leads to: 



6x w 

— * 0.89^- 

Xw 



We now turn to the case of a non-vanishing fa. We can write: 

/iV d 3 Sc 2 



y8 3 



hence for a given resonance we have: 



ft 



3r 



dr 3 



l 



m — - I — Xm + _ X,,, + , X„ 



6r 



(32) 



(33) 



(34) 



from which we infer an approximate expression of the relative shift of the wake: 

±Q.19fah (35) 



6x w 



fah 
~ 3 



-l 1 A 1 / 

u + -A + -uA' 
2 6 
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We can evaluate the impact of those shifts on the differential Lindblad torque using 
the calibration presented in appendixO We use the approximation that: 

xoj-z = -AT<w> (36) 

0*0(7) 

with A » 2, where /"o(/> is the one-sided Lindblad torque exerted at outer (inner) 
Lindblad resonances, and where Xo(r> is the distance of the outer (inner) wake to 
corotation. We eventually write the outer and inner torques respectively as 

r = -r°(i + i/A) (37) 

// = r°(l-J7ft), (38) 

where r° is the one-sided Lindblad torque (average of the absolute value of the in- 
ner and outer torques), and where 77 is a dimensionless parameter that quantifies the 
degree of asymmetry between the inner and outer torque. We have the following ex- 
pression for the one-sided Lindblad torque: 

r° = OAZGfyctfr 3 , (39) 

where the numerical coefficient is easily obtained from numerical simulations. Both 
two dimensional runs with a softening parameter e = Q.6H and three dimensional 
runs with a very small softening parameter yield the same coefficient of order 0.4. 
We make the following comments: 

- the coefficient rj depends on whether the disk considered is two or three dimen- 
sional, and in the former case it also depends on the softening length of the poten- 
tial. In our set of calculations, it turns out to be of order 77 st 2.2, but this value is 
highly specific to our case. As we explain below, nevertheless, this does not have 
any impact on the final result. 

- As P2 varies the outer and inner wakes both recede (or approach) corotation, as 
can be seen from Eq. d3"TV We therefore expect a minute effect of those shifts on 
the differential Lindblad torque. 

- On the contrary, terms arising from yield shifts of outer and inner torques that 
have same sign, thus giving a much stronger (cumulative) effect on the differential 
Lindblad torque. The effect in that case will essentially scale with the one-sided 
Lindblad torque, and the asymmetry term T]h does not feature at lowest order. 



Using Eq. 136K we have: 

8r = -ATo— (40) 

x 

5xj 

sn = -ah— (4i) 

Xi 

We consider separately the case of terms of even order, which have Sxo/xo = Sxj/xr, 
and the case of terms of odd order, for which dxo/xo = —Sxi/xi. In the first case we 
find 

sn = -An—, (42) 

x 
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Fig. 1 Lindblad torque dependence on the second and third derivatives of temperature. Solid lines stand for 
the case of the third derivative (fis). The curve with diamonds stands for a reference calculation (a\ =3/2 
and E = 2H). The curve with triangles has a>i = and E = 2H, whereas the curve with squares corre- 
sponds to ai = 3/2 and E = 3H. The dashed lines show the dependency of the differential Lindblad torque 
on the second derivative of temperature ifii): the line with diamonds represents the results of numerical 
simulations, whereas the plain dashed line represents the expectation of Eq. {48). Some runs with high val- 
ues of |/32 1 or I/S3 1, which have extremely high or extremely low sound speed near the boundaries, happen 
to have crashed in the very first time steps, which explains why some points are missing in this plot. The 
dot-dashed line shows the linear dependence on ,83 of the fully analytic torque estimate, extrapolated to 
larger values of /S3, for the case a\ =3/2 and E = 3H (see main text for details.) 



whereas in the second case we get: 

6 r L = 2Ar () ^. (43) 

x 

In this last case, using Eq. ( 133b . we are left with: 

5r L * 3Air Q /3 3 h (44) 

which can be recast as: 

6r L « 1.2703/^ (45) 

The numerical calculations presented in Fig. Q] exhibit a dependency in broad 
agreement with the estimate of Eq. (l45l l. as a linear regression fit on this data yields: 

5r L = K/3 3 r Tei , (46) 

where the actual value of K, comprised between ~ 0.75 and ~ 1.0, depends on the 
details of the setup, and essentially on the value of E (see Fig. [TJ. This figure also 
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shows the dependency on fa of the analytic torque estimate (see appendix [B] for de- 
tails). We note that our fully analytical estimates do not cover the whole range of fa, 
as in our simple implementation the algorithm that seeks the location of a resonance 
fails for very distorted profiles. Rather, we find that the Lindblad torque value has a 
linear dependence on fa for low values of fa (^1 < 0.3), and that its slope coincides 
quite remarkably with that obtained in numerical simulations (we find a slope K of 
the analytic estimate of 0.71 for E = 2H, and of 0.95 for E = 3H). It is surprising, at 
first, that such a good agreement is obtained between numerics and analytics, whereas 
the analytic estimate performs poorly for power law disk profiles (see section|4]i. We 
note, however, that the further a resonance lies from the orbit, the larger its shift, so 
that most of the torque change is due to resonances of order lower than those which 
predominantly contribute to the torque. This corroborates the fact that the inaccuracy 
of the fully analytic torque expressions likely comes from highest order resonances, 
as a consequence of the approximate form of the torque cut-off. The good agreement 
that we find also justifies, indirectly, our substitution from Eq. ( TT3l > to Eq. ( fT4t . 

It is also noteworthy that this comparison between numerics and analytics com- 
pletely fails if one does not include the factor Q/k in Eq. d9TT >. In that case one finds 
a nearly flat dependence of the torque on fa. 

We now work out the torque change due to fa. Using Eqs. OTb and d42l) . we are 
led to: 

(uAy 1 + l+ uA/2 
5r L = -Ar L — -j3 2 h » -337J3 2 hr L (47) 

Using the fact that, in our case, we turn out to have 2r re f (the actual coefficient 

depends on the qt^'s and /3,'s), we eventually write: 

r L ~ 2.7 WW ~ 0. 14y3 2 r ref . (48) 

The dependency on fa, as expected, is much weaker than the dependency on fa, so 
we shall neglect it. Figure Q] shows that the actual dependence of the Lindblad torque 
on fa is even weaker. 



6.3 Higher order derivative of surface density 

We now turn to the dependency of the differential Lindblad torque on the higher order 
derivatives of the surface density, assuming that the temperature profile remains a 
power law of radius. For this purpose we recast Eq. (l20l as: 



t t -> 1 d log 2" t5 2 1o2 2" 4 I d r c^\ 

k 2 =Q 2 k + (2- fa)Q 2 K h 2 —5- + Q\h 2 — g + -d r c] + rd r ^ . (49) 
dlogr d(logr) 2 r \ r ) 

In a similar fashion, we transform Eq. (fT8l as: 



rU = rL2 v + rU K h — V o r c„ 

K K dlogr 

so we can write: 



(50) 



6D = (2-fa± m)Q 2 K h 2 —2- + Q\h 2 S (5 1) 

ologr d(logr) z 
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Considering as previously that terms that contribute significantly to the torque have 
m » 1, Eq. dTTb leads us to: 

h 2 rd5\ogS h 2 rd 2 6\ogZ 

6r »< = -i — HI ± T~ ~H7i ( 52 ) 

3 ologr 5m a(\ogry 

As in the case of perturbations of the sound speed profile, we entertain separately the 
cases a2 + and a?, + 0. 
In the first case, we have: 

6rm ^( x± L\ (53 ) 
3 V ml 

hence 

S JH = ^ [i + (^)-ij „ o.76/ M2 (54) 

x m J 

This dependency is almost the same as the one predicted by Eq. (fJTJ on /?2, and the 
torque has nearly same dependence on a,i, namely: 

Sr L * 0.12a 2 r ref . (55) 

In addition to the effect of the resonance shifts there is an additional torque variation 
arising from the change of surface density itself. The perturbed surface density corre- 
sponding to a2 being even in x = r - r c , its impact on the differential Lindblad torque 
is very small and is not contemplated here. 

We now evaluate the torque variation due to the third derivative of surface density. 
We firstly work out the resonance shift, and then evaluate the additional effect arising 
from the torque weighting by surface density. 

We write: 

86 log X 1 x 2 ff3 



and 



hence Eq. d52l > gives: 



<91ogr 2 r^ti 1 

d 2 5\ogZ x aj, 
<9(log r) 2 r c h 2 ' 



(56) 



(57) 



2 

X X 



^te^r (58) 



6x m h I _j 1 



±— \ u + —A I Qf3 ~ 0.73/zqt3 (59) 



so that: 

Xm 3 y 2 

A calculation similar to that of Eqs. (PHT i to (l45l l yields: 

sr L * 1.23r ref a 3 (60) 

In addition to the resonance shift, the change of surface density also impacts the 
differential Lindblad torque. One can easily see that the surface density contribution 
tends to compensate the wake shift. Assume that 013 is positive. The resonance shift 
given by Eq. d59l shows that the wake is shifted outward, which yields an increase 
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of the inner torque and an increase of the outer torque (decrease in absolute value), 
hence an increase of the differential Lindblad torque, as seen previously in the case 
of a non-vanishing ^83. However, there is now also a larger surface density at the 
outer wake, and a lower surface density at the inner wake. The contribution of this 
surface density change is opposite to that of the resonance shift, hence we expect the 
dependence of the differential Lindblad torque on a 3 to be weaker than that on (3^,. 
We show hereafter that the contribution of the surface density change is of same order 
than that due to resonance shift (albeit with the opposite sign), and then we resort to 
numerical calculations to assess the actual dependence. 

In order to evaluate separately the impact of the surface density change on the 
one-sided (then differential) Lindblad, we use the torque expression of Eq. (|9Qt , so 
that we can write the relative torque change as: 

fo = y-p . (61) 

Using Eq. ((6]), we recast Eq. doTt as: 



6r° a 3 h Zm x 3 m r, 



6 2 m r n , 



(62) 



where 



r m - a 2/1 
Xm = __, ± _y__ + 1 . (63 ) 

In Eq. d62l . whenever a value of |Jc m | is larger than E/H = 3, the half size of the region 
over which a 3 is constant and non-vanishing, we clamp it to +E/H. This happens for 
the smallest values of m and is meant to adhere strictly to the set up of numerical 
simulations. The evaluation of Eq. d62l yields 

— = U6a 3 h (64) 

We deduce that the surface density impact on the differential Lindblad torque is, using 
Eq. (HI: 

Sr L = 6r + Sn = l.l6a 3 h(r - /» « -232a 3 hr° ^ -0.93a 3 r rei . (65) 

The net effect of the surface density change, which combines the effect of the reso- 
nance shifts, given by Eq. d60l >. and the change of the background density itself, given 
by Eq. d65b . amounts to: 

5r L * 0.3or 3 r ref (66) 

We see therefore that the effect of the third order derivative of the surface density 
should be weaker than that of the third order derivative of the temperature. Numer- 
ical calculations presented in Fig. [2] show that the net dependence is actually much 
weaker, which shows that the two effects described respectively in Eqs. (loTJt and (l65l l 
approximately cancel out. As a consequence, for our purpose we can consider that the 
differential Lindblad torque does not depend on a 3 (as stressed in section I5T1 the an- 
alytic estimate of the shift of resonances that most contribute to the torque is used as 
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Fig. 2 Differential Lindblad torque dependence on (solid lines) compared to the dependence on /?3 
(dashed line, which corresponds to the reference calculation of Fig. [TJ. The differential Lindblad torque 
turns out to be essentially insensitive to a^. Besides, its weak dependence on 0:3 is not a simple linear 
relationship, and its shape depends (slightly) on the value of ct\ (ct\ = for the curve with squares, 
whereas ai = 3/2 for the curve with triangles). The dot-dashed line shows the dependency on 0:3 of the 
full analytic estimate of the torque, extrapolated from lower values of Qa^l < 1, see section l6~2l . 



a guideline and the actual dependency found in numerical calculations is favored for 
the reasons exposed in that section.) Fig.|2]also shows the dependency of the full ana- 
lytic torque estimate (see appendix[B]for details) on (X3. The dependence has opposite 
sign that the one inferred from the simple calculation based on the wake shift. The 
surface density weighting again overshoots the effect of the resonance shifts. This 
observation is reminiscent of the results for power law profiles given in section [4] 
where the coefficient of a\ of the fully analytic torque estimate was found to have a 
sign opposite than that of fi\ . 

We have seen that the only non-negligible contribution to the differential Lindblad 
torque from higher order derivatives of the density or temperature is the term that 
stems from the third order derivative of the temperature. 

We can now propose a tentative generalized differential Lindblad torque expres- 
sion in non power-law disks, for the two dimensional disks with softening length 
e = 0.6//. Out of the four dependencies that we have worked out (on [$2, /?3, &2 and 
a-i), the dependency on given at Eqs. (f43T>-(l46b. is much larger than the other ones, 
given by Eqs. (1481 1. d55l ) and d66t . It is therefore the only contribution that we shall re- 
tain to evaluate the Lindblad torque at locations where the disk profiles are not power 
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laws, so that our tentative generalized differential Lindblad torque reads: 

r L = -(2.oo-o.i6 ffl + i.npx -fc)r mU (67) 

which should be put in contrast with Eq. ( fTOb . 

Eq. d67b raises the question of whether one should push the expansion further so 
as to include higher order derivatives. This question can be circumvented in a manner 
that is exposed in the next section. 



7 Torque behavior at a temperature transition 

7.1 Disk model 

We now consider a disk model in which the surface density is a power law of radius 
(E oc r~ 3 / 2 ), but the sound speed (hence the temperature) exhibits a more complex 
radial dependence given by: 



c s (r) 



1 



(hi + h 2 ) + (hi - h 2 ) tanh 



G K (r)r, 



where w, the width of the transition, is given by: 

w = r t ho 



(68) 



(69) 



and where the disk aspect ratios beyond and prior to the transition are respectively: 

hi = h x 1.1 1/2 (70) 



h 2 



h x i.r 
h x i.r 1 ' 2 



In Eqs. ( f69b and d70b the fiducial aspect ratio ho was chosen to be ho = 0.04. The 
parameters above determine a sound speed jump of ~ 10 % over a scale length ~ 
2hor, corresponding to a temperature jump of ~ 20 %. 

We note that the dust Rosseland opacity k r can exhibit discontinuities of a s much 
as a factor of two as a function of temperature (Henning and Stog nienkol 1996). since 
the mid plane temperature (which matters for torque calculations) scales with T e ff k 1 J 4 , 
where T e a is the e ffective disk te mperature, which has a smooth profile imposed by 
the accretion rate (ICassenl 1 19941) . we can expect jumps of the mid plane tempera- 
ture of as much as 2'^ 4 ~ 1.19. We therefore regard our simplified disk model as a 
bona fide representation of disk properties at an opacity transition, albeit somehow 
extreme. 



7.2 Lindblad torque 

Figure [3] shows the results of numerical simulations, in which we sample the tran- 
sition opacity with 51 different runs. The run number i (i e [0,50]) considers a 
planet of mass M p = 10~ 6 M*, orbiting on a fixed circular orbit of radius r,- = 
[0.4 x (//50) + 0.8]r,. 
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Fig. 3 Lindblad torque near the opacity transition, as a function of the orbital radius. The solid line shows 
the data from the numerical simulations (each diamond stands for a given ran), the dashed line shows the 
results of Eq. HO) and the dash-dotted line shows the results of Eq. {67j. The dotted line shows the result 
ofEq. {79) 



We note on this figure that the Lindblad torque displays a negative peak at the 
opacity transition, whereas the prediction of the standard formula for power law disks 
exhibits a positive peak (as a straightforward consequence of the decrease of fi\ at the 
transition). Finally, Eq. (|67| > does reproduce a negative peak, as expected, but with an 
exaggerated amplitude. 

The origin of the discrepancy can be traced back to the conditions under which 
Eq. d67l i apply. When deriving this equation we assume that the third derivative of the 
temperature profile is constant over the whole region of tidal interaction between the 
planet and the disk, that is at least over [r p - 2H, r p + 2H]. This is not the case of the 
present model, in which the transition occurs on a short spatial scale of order ~ H. 

It is therefore more accurate to keep the local values of the first and second order 
derivatives at the outer and inner wake, in Eq. (|28K which is therefore recast as: 



Sr. 



mh r r 



'sQ. 



(71) 



where c 2 s0 denotes the sound speed at the orbit, and where the superscript + means 
that the value must be taken either at the outer wake (+) or at the inner wake (-), and 
where the c superscript means that the value must be taken at the planetary orbit. We 
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have made use in Eq. (fTTT i of the approximate equalities: 



(72) 



As previously, we leave the exact location at which these values are to be evalu- 
ated as a function of a free dimensionless parameter D, so that variables relative to 
the outer (inner) wake are to be taken at r c + DH (r c — DH). 

Using Eq. ( 136*1 ). we have: 



AT 

Sr L - — - (8xo + 6x1), 
H 



(73) 



where the shift of the outer wake 5xo and that of the inner wake 6xj, using Eq. ( fTTT i. 
are given respectively by: 



8x0 - 
8x1 = 



#2 

3r 
3r 



^ -(fit -ft) 
ft 

^f+^-ft) 



so that the total variation of the Lindblad torque amounts to: 



sr. 



Ar°h 



+ (2ft-ft-ft) 



(74) 
(75) 

(76) 



We note that both terms in the bracket of Eq. (176) stand for finite difference approxi- 
mations of the third derivative of the temperature, as: 



ft -ft* 2D ft 
2ft -ft -pi *D 2 ft 



(77) 
(78) 



Performing this substitution, we would be led to an equation similar to Eq. d45l l. As 
pointed out above, however, we prefer to keep the actual values on either side of the 
orbit, for ft itself may vary significantly over the region of excitation of the wake. 
We note that the coefficient of ft, not surprisingly, is of the order of magnitude of 
the negative of the coefficient of fi\ in the standard torque formula. The last term of 
Eq. (l76"T l therefore amounts to substituting, in the standard torque formula, the value 
of ft by the average of ft and ft. This is rather intuitive, as it is the pressure gradient 
at these locations that shifts the resonances. In a power law disk, these three values are 
equal and the standard formula applies, but if there is a non vanishing third derivative 
of the temperature, {ft + ft) 12 may differ from ft. 

This effect, nonetheless, has only the consequence that it smears slightly the pos- 
itive peak of Fig. [3] Most of the effect is actually due to the first term in the bracket 
of Eq. 1761 . which comes from the second order derivative of the pressure, that con- 
tributes to change the epicyclic frequency and yields most of the resonance shift. 
Since the third derivative of the temperature is not constant over the torqued region, 
the relative shift of the resonances with different values of m differs from the case of 
section I0T2I with the consequence that the coefficient of the dependence on the values 



of (/?2 — ft) may differ, so that we leave it as a free parameter in our final fit. 
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Incorporating these two changes into Eq. ([TUt . we eventually obtain the following 
Lindblad torque expression: 



r, = - 



2.00 -0.16a! + 1.11 1 1 -O.80g^-/3 2 ) 



r ref , (79) 



where we remind that a quantity with a + (-) exponent must be evaluated at r — 
a + DH (r — a - DH), and where D is a free dimensionless coefficient that we 
adjusted to get the best agreement with numerical simulations. We find that we get 
the best agreement for D = 0.2, and in that case the coefficient of - /3 2 must be 
0.8. We remind that @2 is defined in Eq. (7). The results of Eq. d79l are displayed 
in Fig. [3] They show a good agreement with the Lindblad torque inferred from the 
numerical simulations. Some comments are in order about the relevance of the fit 
given at Eq. ( 1791 : 

- In the case in which the third order derivative of the temperature (/J?) vanishes 
over the torqued region, and therefore in particular in power law disks, Eq. ( 1791 
reduces to the standard formula of Eq. (fTOt . since in that case pi - fi 2 = 0, and 
(jS| + p\)/2 = ySj . Eq. ( 1791 can therefore be used everywhere in the disk. 

- As the planet moves radially and samples different locations of the opacity tran- 
sition, the relative shift of the resonances with different values of m varies, so that 
a relationship such as the one proposed can only be approximate. We therefore 
seek agreement where the effect is strongest, that is at the nominal position r — r, 
of the transition, where the peak of the Lindblad torque is observed. 

- We shall see in section l731 that the expression proposed satisfactorily accounts for 
the torque measured for transitions of different widths, and that the strength of the 
effect scales the inverse cube of the width. Transitions significantly more narrow 
than those contemplated here are Rayleigh unstable (lYang and Menoul 120101) . 
whereas the effect will be negligible in transitions significantly wider. Eq. (T79b 
constitutes therefore an acceptable description of the effect for all situations in 
which it is measurable. 

- Contrary to what has been done with the temperature slopes, the surface density 
slope coefficient is left unchanged. There are several reasons for this: (i) the coef- 
ficient of a\ is small, hence no measurable change is expected if one replaces its 
value with the average of its values at the outer and inner wake; (ii) we have seen 
in section [63] that the additional effect of surface density weighting of the torque 
tends to cancel out the effect of the resonance shifts, hence we do not write a term 
in at - a 2 , similar to that in (3^, since that term would have a small coefficient, 
and (iii) at an opacity transition, no significant effect is expected for the surface 
density, counter to the temperature. 

- The changes with respect to the standard formula consists in: (i) an evaluation 
of the shift of the resonances using the slope at the wake location itself, rather 
than at the orbit (this corresponds to the term in fi\ + fT x and it is a small effect), 
and (ii) an evaluation of the change of the epicyclic frequency itself, the value of 
which defines the location of Lindblad resonances (this corresponds to the term 
in ySj -jSj, and it is the most important effect). From Eq. ( l20l . we see that the 
derivative of the temperature to second order suffices to evaluate the change in the 
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Fig. 4 Difference between the true value of the torque and the one inferred from the standard torque 
formula, as a function of the distance to the nominal radius of the transition. The solid line represents the 
results for numerical simulations (each diamond represents a run), and the dashed line shows the results 
for the analytical estimate. The inset plot shows the total Lindblad torque as a function of the distance 
of the nominal radius of the transition, obtained using a fully analytic estimate (dotted curve) or using 
the corresponding standard torque expression given by Eq. (1121 . We recover the behavior found from 
numerical simulations at Fig. [3] namely that the actual profile is wider, and displays a reversed peak at the 
top. The dashed line of the main plot is essentially the difference between these two curves. 

epicyclic frequency, and no further expansion is needed, provided the derivatives 
are evaluated at a location intermediate between the orbit and the wake, rather 
than at the orbit. 

7.3 Comparison to the fully analytic torque expression 

We show in Fig. [4] the difference between the actual torque value and the value in- 
ferred from standard torque estimates, both numerical simulations and for fully ana- 
lytic estimates. Differently said, this differences assesses the magnitude of the effect 
of higher order terms of the temperature profile. Although the global shape of the 
curve of the fully analytic expression is broadly similar to that obtained from simula- 
tions, it displays a larger amplitude, which is consistent with the fact that the standard 
torque expression obtained by summing torque estimates on all resonances, given by 
Eq. ( fl2l i. has a larger coefficient for fi\ than the standard torque expression given by 
numerical calculations, given at Eq. (TTOt - We had noted in section [6721 the excellent 
agreement between the linear estimate and numerical simulations for the dependence 
of the torque on /?3. Our tentative explanation was that the effect was primarily due to 
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relatively low order resonances, which had larger shifts. Here, the resonances that are 
most shifted are the ones which lie the closest to the orbit, i.e. those of higher order. 



7.4 Corotation torque at the transition 

We now contemplate the other torque component, namely the corotation torque, 
which has been thus far omitted. We examine here the vortensity related corotation 
torque, in the linear limit. We therefore focus on the indirect effects of the tempera- 
ture profile: it alters the rotation profile and therefore the vorticity profile, hence the 
vortensity gradient. 

We wish to evaluate the dimensionless coefficient 

^ = rflog(^) 
dlogr 

Since the profile of Z is a power law, we focus on the part related to the vorticity 
a>, given by Eq. (0. We use Eq. and we keep only the derivatives of the sound 
speed, which varies over the length scale H, whereas other quantities vary over r. 
This yields: 

d\c] 

to—is. (81) 
Considering that a> ~ Q/2 everywhere in the disk, we are led to: 

*V — fa- (82) 

Eq. ( f82b is an approximate relationship that holds wherever there is a significant peak 
of /?3 in the disk. In the power law parts of the disk, we have *V = 3/2 — a\. We can 
deduce the corotation torque at the transition using Eqs. (TT4l and (l82l . This gives: 

r c ~ o.6(fii -j3 3 )r ref (83) 

Notwithstanding its standard dependency in /Si , the corotation torque should therefore 
display a boost opposite to that of the Lindblad torque at an opacity transition, as we 
can see by comparing Eqs. ( l83l and (l46l >. Although it may seem at first glance that 
the balance is in favor of the Lindblad torque, such is not the case for the following 
reasons: 

- the previous section has shown that the peak of Lindblad torque observed at the 
transition was less than that predicted by Eq. (|46*1 i. This was understood as being 
due to the fact that the third derivative is not uniform over the torqued region. 

- In contrast, it is the third derivative of the temperature at the orbit that gives the 
corotation torque. 

- In addition to this, we remind that the torque value in these two dimensional cal- 
culations is very smoothing dependent. The value of the smoothing adopted here 
(0.6//) is adapted to the Lindblad torque, in the sense that it provid es a one-sided 
Lindblad torque of t he order of that in a three-dimensional disk dMivoshi et all 
1999; iMassen , 120021) . However, this value of the smoothing tends to underesti- 



mate the corotation torque, for which a smoothing length of ~ 0.4// should rather 
be used (Casoli & Masset, in prep.) 
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Fig. 5 Total torque near the opacity transition, as a function of the orbital radius. Black lines show the 
results for the standard transition whereas grey lines show the results for the narrow transition. See the 
main text for the signification of each line style. 



- Finally, the corotation torque, depending on the planet mass, on the disk thickness 
and on the amount of viscous diffusion, eventually remains linear at all times, 
or settles into the so-called horseshoe regime. The horseshoe drag bears same 
dependencie s on the disk parameters as the linear corotation torque, but it has a 
larger value (Paarde kooper and Papaloiz ou 2009). Only for very small values of 
the disk's effective viscosity does the horseshoe drag enter a so-called saturated 
regime, so that the total torque e ssentially amounts to the differential Lindblad 
torque ( Masset and Casoli , 2010l) . 



7.5 Total torque at the transition 

We present the total torque at the opacity transition in Fig. [5] Results are presented 
for two transitions: the same as the one contemplated in previous paragraphs, de- 
fined by Eqs. d68l > to (l70l . which has a width r,ho ~ H, and a more narrow one (w f = 
Q.lr t h .Q ~ 0.7//), which is nevertheless still largely Rayleigh stable dYang and Menoul 
2010h . The solid curves with diamonds show the results directly obtained from nu- 
merical simulations (each diamond stands for a different run). The dotted lines show 
the expectations of the standard total linear torque expression of Eq. (fTTI i. and the 
dot-dashed lines show the sum of Eqs. (TT~0t and (TBI , that is the standard total torque 
expression, except that we use the proper value of 'V instead of 3/2 - a\ . Finally, the 
dashed lines show the sum of Eqs. (fl4l i and (l79l . 
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We see that the standard torque formula fails at reproducing the strong peak of 
the total torque. Even in the narrow case, it does not predict the observed torque 
reversal. On the contrary, the mere addition of Eqs. (fTOl and (TBI overshoots the 
real torque, essentially because it disregards the Lindblad torque sensitivity on higher 
order derivatives of the temperature. Finally, the dashed lines, that corresponds to the 
sum of the proper expressions of the Lindblad and corotation torque, are in excellent 
agreement with the torque measurements. They correctly predict the occurrence of 
torque reversal, and the magnitude of the peak at the transition. 

We also note the steep dependence of the peak on the transition width w t . The 
peak value scales indeed with y03, which is proportional to w^ 3 . 



8 Discussion 

8. 1 Extension of results to more general cases 

We have seen that both the Lindblad torque and the corotation torque have a depen- 
dency on the third derivative of the temperature, with opposite signs. The balance is 
in favor of the corotation torque, which, even in the linear regime, could be enough 
to stall migration at opacity transition, provided the transition is sufficiently narrow 
(we comment that our calculations, which exhibit such stalling, are conservative in 
the sense that they underestimate the corotation torque by the use of a large softening 



of the potential). Early work on this kind of traps (see Masset etal 20061 at section 
5) did not evidence this possibility of stalling, essentially because the transition was 
too large (~ AH) and the potential softening value even larger (e = 0.7//). 

In the general case, one should take proper care of the value of the corotation 
torque by incorporating possible non-linear effects, namely the horseshoe drag regime 
and its partially saturated state. Details about this kind of procedure can be found in 



Masset and Casolild2010h 



A result similar to the one presented above is expected to hold for three dimen- 
sional disks. Namely, the two effects exhibited in section Q should be introduced in 
the same manner as for two dimensional disks: 

- The temperature gradient at the orbit should be replaced by its average at the outer 
and inner wakes. 

- A term in ySj - fi^ should be added. We remarked in section 16.21 that the one 
sided Lindblad torque has same value for three dimensional disks and for the two 
dimensional disks that we considered. Since the effect arising from the change of 
epicyclic frequency is a cumulative effect (it has same sign for the inner and outer 
torques) that scales with the one-sided Lindblad torque, we expect the coefficient 
of this term to be the same as in our two dimensional disks, that is to say ~ 0.8. 

Still, a standard Lindblad tor que formula i n thre e dimensional disks has not been 



worked out yet. The work of Tan aka et all (120021) is limited exclusively to globally 



isothermal disks. The calibration exposed in appendixlAlis hardly achievable in three 
dimensional calculations. In particular, one cannot resort to nested meshes centered 
on the planet, as the low resolution in the parts of the horseshoe region that lie far 
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from the planet would act as a source of diffusion, which would impede a full sat- 
uration of the cor o tation torque. The tentative Lindblad torque formula given by 



Masse t and Cas oli (2010) at Eq. (160) was based upon three dimensional calcula- 
tions in the linear regime for which the total torque exhibited a dependency in 0.4/3i, 
and the same dependency was assumed for the Lindblad torque by ignoring a possi- 
ble dependency of the linear corotation torque on f}\. The present work stresses that 
the linear corotation torque has a dependency far from negligible on the temperature 
gradient. Work is in progress to assess this dependence in three dimensional disks 
(Casoli & Masset, in prep. ). It is likely that th e Lindblad torque has also in that case 
a steep dependence on f}\ ( Paardekooper et al , 2010h . 

We stress that our additional term is given in order of magnitude only, but it 
is straightforward to apply, and it yields a value quite close to that obtained from 
numerics in the regions of the disk where the profiles cannot be approximated by 
power laws. 



8.2 Comparison with the analysis of lMenou and Goodman! (2004) 



In their analysis of the migration of low-mass protoplanets near opacity transitions, 
iMenou and Goodman! |2004) found a very strong slowing down of migration, cor- 
responding almost to a cancellation of the Lindblad torque at some locations in the 
disk. On the contrary, we find here a relatively mild effect on the Lindblad torque, 
even though we have considered a rather extreme opacity transition, both in terms of 
narrowness and temperatur e jump, and even though we c onsider a smaller softening 
length of the potential than lMenou and Goodman! (I2004r) . In order to understand the 



reason for this discrepancy, we have tried to reproduce the results o flMenou and Goodman 
(2004) 'or power law profiles and with a sof tening length of the potential equal to 
H. These results were subsequently recast by iD'Angelo and Lubow (1201 Oh . in their 
Eq. (8), as: 



r 



<6=H) 



-(0.80 - 0.77ari + 1.12jSi)r ref . 



(84) 



The most striking feature of this expression is the small value of its constant coef- 
ficient (0.80). This makes the torque value prone to cancellation even for moderate 
values of a\ or j3\. Using our analytic approach detailed in appendix IB] we obtain 

r { lJ H) = -(2.53 + 0.32a! + 0.83/3, )r ref , (85) 

w hich differs substantially from Eq. (l84l . There are some differences in the approach 
of Me nou and Goodm an (2004) which may account for this difference: 

- They use a softening length equal to the disk thickness at the location of the 
resonance, rather than at the location of the planet. 

- They use the adiabatic sound speed in all terms related to non-axisymmetric dis- 
turbances (it is however unclear whether they do use it: th eir figure 2 shows a 
virtually perfect agreement with the results of Ward (1997), when they normal- 
ize the torque to the same quantity, related to the isothermal sound speed - their 
equations (15) and (20) -; one would expect the torque value to be reduced by a 
factor of y — the ratio of the specific heats — if the adiabatic sound speed were 
used). 
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If we incorporate these changes into our analytic approach, we are left with the fol- 
lowing torque expression: 

r Li H) = ~ (lM + °- 24a i + 0.83/?i)r ref . (86) 

We have tried other combinations (e.g., using a radially varying smoothing length 
and using the isothermal sound speed), but none of them has yielded results close to 
Eq. ([84-b . In all cases the constant term remains of order ~ 2. The reasons for this 
discrepancy are unknown. Their approach differs in that they use a continuous torque 
density, which may have a lower accuracy for the relatively thick disk (h = 0.07) 
that they consider, and for the large value of the potential softening length, as only 
relatively low order resonances contribute to the torque in that case. In the same 
vein, the transformation of the Laplace coefficients to modified Bessel functions of 
the second kind may further degrade the accuracy, but none of these explanation 
sounds plausible at accounting for the large difference between Eq. ( l84l l and Eq. ( l86l l. 
We believe that the huge effects that they get at opacity transitions are linked to the 
smallness of the constant coefficient in Eq. ( l84t . Our results, both from numerical 
simulations and from analytic estimates, are quite different: we find that the effects 
of opacity transitions (even extreme ones) on the Lindblad torque cannot reduce it 
by more than a factor of about ~ 2. However, the total torque can still be reverted 
at these locations, not by the Lindblad torque but by the corotation torque that can 
undergo a boost originating from the strong departure from Keplerianity of the flow 
at the transition. 



9 Conclusions 

We have obtained an expression for Lindblad torques, given at Eq. d79l l, that general- 
izes the standard expression of Eq. (TT~0b wherever the disk profiles are not power laws, 
in particular at opacity transitions. Essentially one additional term plays a role, at lo- 
cations where the third radial derivative of temperature does not vanish. The main 
effect comes from the alteration of the epicyclic frequency, which depends on the 
flow's vorticity. The latter differs significantly from its nearly Keplerian value wher- 
ever the pressure gradient varies over a short length scale. The resulting variation of 
the epicyclic frequency shifts the planetary wake and yields a substantial change of 
the Lindblad torque. 

Our expression has been obtained for locally isothermal disks. In disks for which 
an adiabatic expression is better suited for the Lindblad torque, the value found here 
can simply be divided by y dBaruteau and Massetil2008l:|Paardekooper and Papaloizou 
120081 : iMasset and Casolii l2010i : IPaardekooper et all 120101) . whereas the value of D 
should be multiplied by y 1 ^ 2 . The transition from the isothermal valu e to the smaller 
adiab atic value, as a function o f the thermal diffu sivity, is described by Masset and Casoli 
(|2010h at Eqs. ( 1 54U 1 561. bv lLvra et all d20 1 Oh at Eq. (10), or by Paard ekooper et all 
( 201 lh at Eqs. (45)-(47). These three procedures are in broad agreement. 

There are side results to this analysis. Firstly, we find evidence from numerical 
simulations that the linear corotation torque in a locally isothermal disk includes a 
term that scales with the temperature gradient, and that the dependence is as steep as 
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the dependence on the vortensity gradient, which is an important fact perse but which 
has received little attention so far. Secondly, we find that the c orotation torque, ex- 
plicitly neglected in the analysis of Me nou and Goodm an (2004), displays a boost at 
an opacity transition, because there the vortensity gradient significantly differs from 
its Keplerian value 3/2 — ac\. This boost may be sufficient to halt migration at opac- 
ity transiti ons. An extremum of vortensity is n onetheless unstable to a Rossby wave 
instability dLovelace et a 

1E223 

iLi et alLl2000l) . and viscous diffusion, either laminar 



or turbulent, acts at spreading radially narrow vortensity features. As this was not our 
primary goal, these effects and their interplay with migration have been disregarded 
in this work. 

The present analysis also disregards any adiabatic effects on the corotation torque. 
Should the flow behave nearly adiabatically on the time scale of horse shoe U-turns, 
the additional entropy related corotation torque should be taken into account. 

Acknowledgements I thank Steven Lubow for pointing me out that the linear corotation torque could 
depend on the temperature gradient. 



A Calibration of the Lindblad and linear corotation torques in 
two-dimensional disks with a softened potential 

For the Lindblad torque calibration we perform long term inviscid calculations in which we let the horse- 
shoe drag saturate, so that the asymptotic torque value is considered as a good approximation to the linear 
Lindblad torque (see section [3}- In order to obtain the coefficients Icq to £2 of Eq. {8) for our disks, and 
similar coefficient k' Q to k\ for the total torque we have ran four different calculations, in power law disks. 
Three calculations should be sufficient, but we have ran a fourth one as a control to get some confidence 
that the Lindblad torque can indeed take the form given by Eq. (8). 

Our calculations were carried out for a disk with h = 0.04 and a planetary mass M p = 10~ 6 M„. The 
mesh extends from R mm = 0.5a to R mlllL = 1.6a, while its resolution is N m d = 1650 and Ng = 1000. This 
resolution is meant to give the same cell size as in the calculations of section[6] but with mesh boundaries 
which lie further from the orbit. The horseshoe zone is therefore covered by approximately 9 cells radially. 
The dimensionless parameter that controls the non- linearity of the flow is M p j(h 3 M,) a 0.015. 

In addition, a time average from t = 2 to 4 orbital periods only of the total torque gives an estimate 
of the linear torque value in each case. We comment that over this time frame the torque can still be 
considered as linear. Indeed, the characteristic timescale for the onset of non-linear effects is t he horseshoe 
libration time, multiplied by twice the disk aspect ratio ( Paardekooper and Papaloizou, 2009), which can 
be cast as: 

Tnl ~ ^(q/h 3 )-" 2 , (87) 

a time s cale which can be identified wi t h the time it takes to execute a horseshoe U-turn ( Baruteau and Masset, 
120081) . Paardekoop er and Papaloizoul feO09h mention that the linear results are valid at very early time 
(about two orbits) for their fiducial calculation with q/h 3 =0.1. For our value of q/h 3 , it should be valid 
over approximately up to 5 orbits, while the establishment of the horse shoe drag regime takes about 30 or- 
bits in our case. The results of one of the runs are represented in Fig. [6] in which one sees the three stages 
of the torque : (i) linear stage, at very early times, up to a few orbits, (ii) unsaturated horseshoe drag, 
reached after around 30 orbits, and maintained up to about 120 orbits, and (iii) the saturation phase, in 
which the torque oscillates about the linear Lindblad torque, which constitutes its asymptotic value at time 
larger than 10 3 orbits. 

The torque values, normalized to the reference value of Eq. (|9j, are presented in Tab. [T] This data 
yields an over-constrained system of six unknowns (yo, y\, yz, y\ r y\, y' 2 ) and eight relations. We add a 
constraint, namely ko = (3/2)(i f [ ~k\ )+k' , to ensure that the corotation torque scales exactly with (3/2)-ai 
in power law disks (while we leave its dependence on the temperature gradient as a free parameter), as we 
need to recast it in terms of the vortensity gradient, which differs substantially from (3/2)— o-i when higher 
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Fig. 6 Total torque as a function of time. The torque is shown with two different scales as a function of 
time: over a long time scale to exhibit its asymptotic value repre sented by a dashed line (black curve, lower 
jt-axis), and over a short time scale to expose the transition (see Paardekooper and Paoaloizou, 2009) from 
the linear corotation torque to the non-linear, horseshoe drag regime (grey curve, upper jc-axis). The light 
grey stripe on the left corresponds to the temporal window 2-4 orbits, over which the torque is averaged, 
and which yields an estimate of the Lindblad torque plus linear corotation torque, represented by the dotted 
line. 





P\ 




r£T/r ref 


3/2 





-1.74 


-1.78 


1/2 





-1.91 


-1.32 


3/2 


1 


-2.86 


-2.30 


1/2 


1 


-3.03 


-1.84 



Table 1 Normalized Lindblad torque (at larger time) and total linear torque (averaged between 2 and 
4 orbital periods) for the four calibration runs. 



order derivatives of the surface density or temp erature do not van ish. We note that the linear corotation 
torque, obtained from the semi-analytic work of Tanaka et al (2002), scales indeed exactly as (3/2) — ai in 
power law, two dimensional disksQ From this data we infer the linear regressions given in Eqs. (10) and 
(TTJ. We note that if we do not impose the constraint given above, we are led respectively to: 

r L = -(2.00 - 0.16a! + 1.12/3i)r rcf (88) 

and 

r tot = -(1.09 + ().46ffi + 0.52/3i)r rcf (89) 



4 Although this dependence is not given explicitly for two-dimensional disks by Tanaka et al 12002), 
one can infer it by subtracting the Lindblad torque from the total torque. 
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in stead of Eqs. UP) and Jilt . The comparison gives an idea of the fitting error, as does the first line of 
Tab. [T] Since there is no vortensity gradient and the disk is globally isothermal, both torque values in this 
row should coincide. Sources of inaccuracies may include: 

1 . the narrowness of the time frame over which we average the torque to get the total linear torque (but 
taking a larger time frame would raise the issue of the corotation torque beginning to be non-linear in 
the presence of a vortensity or temperature gradient); 

2. numerical diffusion of the vortensity within the horseshoe region, which may act as a physical viscos- 
ity and may leave a tiny corotation torque at larger time (almost fully saturated, but not quite); 

3. Partial reflection of the wake on the mesh boundaries, which is present at larger time, but not at earlier 
time when the linear total torque is evaluated. As non-reflecting boundary conditions are implemented 
in FARGO, the amount of reflection should be small. 



B Fully analytic torque estimate 

Our method for estimating the Lindblad torque for each azimuthal component of the planetary forcing po- 
tential is essentially based upon that of Ward ( 1997), with two amendments described below. The Lindblad 
torque at a Lindblad resonance of mth order is given by 



"' r m (dDldr) rm ' 
where r„, and D(r) are defined respectively at Eqs. \\5\ and | |16K and where 



(90) 



where 



™ _ rd(p m /dr + 2(Q/K)mf<f> m 

with i; = mc s (r)l rK(r), where t he epicyclic freq uency *: is defined at Eq. U9K and where / = m[Q(r) - 
Q p ]/Q(r). Unlike [Ward (1997), and following Menou and Goodman 12004), we have incorporated to 
Eq. j9U the ratio Q/k. The latter turns out to be unimportant in disks with power law profiles, but it is 
crucial to include it in disks with abrupt changes, where the epicyclic frequency may significantly differ 
from its nearly Keplerian value. Our second modification with respect to Ward's formulation is that we in- 
clude the potential softening length in the evaluation of the amplitude of the mth order Fourier component 
of the disturbing function: 

GM„ ,„ 

*m = -^b'D^ir/a), (92) 
2 C cos m8d6 

b 7l2MaMl ci )=- =■ (93) 

l/AWa) „ J Q _ 2(r/a) Q ~ (r/a)2 ~ (f/a)2 

The form of Eq. 1911 . together with the definition of r,„ and the definition of the forcing function of 
Eq. J93K determines how the torque behaves in the limit of high order resonances, and therefore incor- 
porates the torque cut-off. In the particular case of power law disks, using our fiducial smoothing length 
e = 0.6//, we find r,„ oc exp(-1.8/im) for larger m. In Eq. (93), the value of the softening length e 
is fixed, and it is thus the same for all Lindbla d resonances . This is consistent with the formulation of 
Paardekooner and Papaloizou (2009) and Paardekooper et al (2010) (see in particular Eq. (3) of the for- 
mer). This also corresponds to the standard implementation of the FARGO code, in which the material 
in each zone of the computational domain is acted upon by the planetary potential, which has a given 
smoothing length, and reacts back on the latter through a potential which has the same smoothing length, 
so as to approximately ensure the f ulfillm ent of the action and reaction law. This is in contrast with the for- 
mulation of Menou and Goodman 1 2004), who adopt a smoothing length that vary with r (more precisely 
their smoothing length is the vertical scale height of the disk at the location of the resonance). 

We further comment that we do not reduce the expression of the Laplace coefficients (modified by 
the presence of the smoothing) to expressions involving modified Bessel functions of the second kind, as 
is a standard practice. We rather evaluate directly the integral of Eq. 193) at each resonance, as well as its 
radial derivative. Also, we do not use the convenient device consisting in working out a torque density, 
which converts a discrete problem into a continuous one (e.g. Ward, 1997). Rather, we evaluate dD(r)/dr 
at each resonance, the location of which is formerly determined using a Newton-Ralphson method. 
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Fig. 7 Total inner (squares) and outer (diamonds) torques (in absolute value) and differential Lindblad 
torque (triangles), for a vanishing smoothing length. The vertical scale is the same for both plots. On 
the left y-axis, the torque is normalized to Tq = nq^Za S^Lh ', whereas on the right y-axis it is normal- 
ized to the reference value given by Eq. l|9j. Both values differ by a factor n/h. The latter normalization 
is to be preferred, for it renders the differential Lindblad torque indepe ndent of the disk aspect ratio, if 
the disk is sufficiently thin. T hese plots should be compared to those of Ward 1 1997) (fig. 3) or those of 
Menou and Goodman (2004) (fig. 2). The dotted line shows the differential Lindblad torque for a smooth- 
ing length e = 0.6//. 



In Fig.^we present our results for disks with power law profiles already studied in the literature in 
order to assess the correctness of our implementation. We also display on these graphs the differential 
Lindblad torque that we obtain with our fiducial smoothing length. 



C Calibration of the one-sided torque dependence on the wake's distance 

We evaluate here how the one-sided torque varies as the wake recedes or approaches the planet. We do 
this by means of numerical simulations in which we maintain the planet on a fixed, uniform and circular 
orbit, but we vary its orbital frequency from ran to ran. Doing so shifts the wake in a known amount. We 
measure independently the inner and outer Lindblad torques, and we translate the planet's frequency offset 
SQp into a wake shift Sx w = —(2/3)aSQ„/Q„. The results are presented in Fig. [8] Not surprisingly, the 
dependence is essentially similar for inner and outer torques. We infer from these calculations that: 

£-2**. (94) 

r x w 

from the slopes measure d in 5x w = 0. The distance to the wake x w is assumed to be H (see e.g. the torque 
density distribution in D'Angelo and Lubow, 2010). 
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